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Abstract 

Observed since the 1950's, the seasonal effect on underground muons is a well studied 
phenomenon. The interaction height of incident cosmic rays changes as the temperature of 
the atmosphere changes, which affects the production height of mesons (mostly pions and 
kaons). The decay of these mesons produces muons that can be detected underground. The 
production of muons is dominated by pion decay, and previous work did not include the 
effect of kaons. In this work, the methods of Barrett and MACRO are extended to include 
the effect of kaons. These efforts give rise to a new method to measure the atmospheric 
K/vr ratio at energies beyond the reach of current fixed target experiments. These methods 
were applied to data from the MINOS far detector. A method is developed for making 
these measurements at other underground detectors, including OPERA, Super-K, IceCube, 
Baksan and the MINOS near detector. 
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1 Introduction 



When cosmic rays interact in the stratosphere, mesons are produced in the primary hadronic 
shower. These mesons either interact again and produce lower energy hadronic cascades, or de- 
cay into high energy muons which can penetrate to detectors deep underground. The temperature 
of the stratosphere remains nearly constant, only changing slowly over longer timescales such 
as seasons (with the exception of the occasional Sudden Stratospheric Warming events observed 
during wintertime at high latitudes [HI). An increase in temperature of the stratosphere causes a 
decrease in density, reducing the chance of meson interaction, resulting in a larger fraction decay- 
ing to produce muons. This results in a higher muon rate observed deep underground [I2I3I4I5I6II . 
The effect increases as higher energy muons are sampled, because higher energy mesons with in- 
creased lifetimes (due to time dilation) are involved. This effect permits the measurement of the 
atmospheric charged kaon/pion production ratio. The rate of low energy muons at the surface of 
the earth is also affected by the temperature because the varying production altitude changes the 
chances of the muon decaying before reaching earth, but this effect is not relevant for detectors 
deeper than 50mwe[|3J (meters water equivalent). 



2 Muon Intensity Underground 

The intensity of muons underground is directly related to the production of mesons in the strato- 
sphere by hadronic interactions between cosmic rays and the nuclei of air molecules. It is assumed 
that meson production falls off exponentially as e^^/^^ where Ajv is the absorption mean free path 
of the cosmic rays and X is the slant depth of atmospheric material traversed. It is also assumed 
that the mesons retain the same direction as their progenitors, that the cosmic ray sky is isotropic 
in solid angle at the top of the atmosphere, and ionization is neglected. These assumptions are par- 
ticularly valid for the large energies of the mesons that produce muons seen in deep underground 
detectors such as Baksan [7J, Super-K [8], IceCube [9] and MINOS far detector (MINOS FD) and 
near detector (MINOS ND) [10]. In this approximation, is constant. Two meson absorption 
processes will be considered: further hadronic interactions, dX/ Km, where dX is the amount of 
atmosphere traversed, and M /iz/^ decay. The fractional loss of mesons by decay is given by 



where p = air density, tq = mean M lifetime (at rest) [^, and Pr is the meson rest frame momen- 
tum. M is either a vr or K meson (charm and heavier meson production doesn't become important 
until ~ 10^ TeV). For an isothermal, exponentially vanishing atmosphere, the atmospheric scale 
height H(T) = RT/Mg. The density p is then related to X by p = X cos 9 /H(T) . The critical energy 
eM, the energy that separates the atmospheric interaction and decay regimes, is given by 
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Since most interactions take place in the first few interaction lengths [fTTI . and to first order 
H{T) ^ Ho = 6.5 km, = 0.115 TeV, ej^ = 0.850 TeV, the differential meson intensity 
M{E,X, cos 6) can be written as a function of X [I2I11II : 

for relativistic M, where Nq{E) is the differential M production spectrum which has the form 
E^'^~^^\ Atv is the nucleon interaction length, and Z^m is the spectrum- weighted inclusive cross 
section moment. This differential equation is straightforward to solve using an integrating factor 
and rewriting [|2I11| : 
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Now that an expression for the production and propagation of mesons through the atmosphere has 
been found, a function describing the production of muons must be found. Muons are produced 
from mesons via the two body decay process M — )• fxu. The rest frame momentum for this decay 
is pr = (1 — ^li/^\i)^M /2 (since the neutrino has negligible mass). The differential flux per unit 
cross section is proportional to the differential flux per energy, which can be written 

dn _ BuiM 

dE ~ 2prPL 

where B is the branching ratio and Pl is the momentum of the decaying particle in the lab frame. 
The muon production spectrum for meson M parents is given by Gaisser [,11,1 : 

P,(E,X,cose)= E '^''^ff^ M{E\X,cose)dE'. (6) 

mesons ^min (Xl-j LJ COS U 

Inserting Eq. [5]into the muon production spectrum (Eq. (6]) gives 

-v ^ m [Ejr,, dE M{E,X,cose) 
V^{E,X,cose) = 2^ — — -/ — , (7) 

X cosOn — rM> JEf, E E 



mesons 



where tm = "^^j/^if • Muons are sampled by detectors at one particular depth, so the production 
spectrum must be integrated over the whole atmosphere to find the energy spectrum of interest. 
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The relevant energy spectrum is written: 



^ = r VJE, X)dX ^ Co X ^-(^+1) { ^ -— + 0.635 , , ^ , 

(8) 

where 7 = 1.7 is the muons spectral index lfT2ll . the branching ratio B(K — )■ i/^/i) = 0.635, and 
B(n —7- z/^/i) ~ = 1. The parameters A^r^A- are constants involving the amount of inclusive meson 
production in the forward fragmentation region, the masses of the mesons and muons, and the 
muon spectral index [fTTIl : 



17+1 



^7r(X) = JZ T — • (y) 

(1 - r^{K)) 7 + 1 

The integral of the production spectrum can be written in the form [|2l 

UE) = / dE,-^. (10) 

This is the total number of muons with energy greater than the minimum required to reach an 
underground detector. The threshold surface energy required for a muon to survive to slant depth 
d{0^ 0) (mwe) increases exponentially as a function of d and parameters a(E) and b(E) IfTTI . Since 
a and b depend on energy, an iterative procedure can be used to find the threshold energy [jT2|: 



Eth = E:^\e, 0) = [e- + ^) e^'^(^'^) - ^, (11) 

where the energy-dependent parameters a = 0.00195 + 1.09 x 10^^1n(_E') GeV cm^/g and b = 
1.381 X IQ-^ + 3.96 X 10"*^ \n{E) GeV cm^/g at column depth die, 0). The threshold en- 
ergy at the minimum depth of the detectors considered in Sec. [T] are shown in Table [Tl Eq. [TOlis 
approximated [2J as 

/ ~C xE-"'( - ^ (12) 

*^ l,7 + (7 + l)l-l^thCOS^/e^^7 + (7 + l)l-l^thCOS^/ex;' 

where 0.635 •yli^/A^-r(K/7r) = 0.054 HIl, lEl and r(K/7r) is the atmospheric K/tt ratio. 



3 Temperature Effect on Muon Intensity 



The temperature changes that occur in the atmosphere are not uniform, instead occurring at mul- 
tiple levels, and neither muon nor meson production occurs at one particular level (see Figs. [H [21). 
The perturbations that variations in temperature cause in muon intensity are small, and as a result 
properly chosen atmospheric weights can be used to approximate the effective temperature of the 
atmosphere as a whole, T^s. Define r]{X) = (T(X) — Tcs)/Tcs, and ejv/ = 6^/(1 + i]), where e^^ 
is the constant value of eM when T = Tes- This is the temperature that would cause the observed 
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Table 1 

Threshold muon energy for current underground detectors. 



Detector 


Min. Depth (MWE) 


Eth(GeV) 


MINOS ND d 


225 


51 


Baksan El 


850 


234 


IceCube B 


1450 


466 


MINOS FD in 


2100 


730 


Super-K H 


2700 


1196 


OPERA aa 


3400 


1833 



muon intensity if the atmosphere were isothermal. To quantify the temperature effect on intensity, 
the temperature dependence of Eq.|2] needs to be considered. The meson production term in Eq.[3] 
(which applies to any charged meson: K, vr, etc) can then be expanded: 
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The analytic solution to this differential equation is difficult to find since r]{X') is an arbitrary 
function of X'. A solution to first order in i]{X') can be found by expanding the exponential in a 
power series, and then following the procedure outlined above, beginning with Eq.[5l This solution 
can be written as X, cos 6*) = + M^,v/heTe M^{E, X) is the solution where cm = e°j, 

which occurs at temperature T = Teg and M^{E, X) is given by: 



A1'(B,.Y,«) = %iiV„(£)e-VA.'(iL) ^ dX'!!^ ' ^ ' 

An \AmJ EcosOJo 



[el^/Ecose + l e%/Ecos9 + 2 2\ elj/E cos9 + 3 "' ^ 




If EcosO ^ then the integrand is very small and ri{X') = ri{X). This is the case when in- 
teractions dominate, as time dilation effects allow these very high energy mesons to travel great 
distances before decaying. If E cos 9 <^ e°/, then the mesons will not travel very far before decay- 
ing and the integrand is large only when X' is near X, and again, 77 (X') can be taken out of the 
integral 

Writing the solution of Ai where T = Teg as and letting = e°j(l + 77), an expression for 
the change in muon production induced by temperature variations can be found. Define AA^ = 
M. — M^, then to first order in r] 
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Using Eq.|7]and Eq.[8l an expression for the change in differential muon intensity can be found: 



where 
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Now, a solution to this integral can be found for E^ 3> t\.i{I^) and for E^ <ti t\i{I^): 
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These expressions can be combined in a form that is valid for all energies (Eq.[8]): 
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The exact solution for I^,^{E^) has been replaced with an approximation that preserves the physical 
behavior of the system at low energies. These low energy mesons are relatively insensitive to 
changes in temperature because they decay before they have a chance to interact. So, this equation 
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describes the expected behavior that mesons at very low energies will decay fairly high in the 
atmosphere. These mesons will not contribute any muons to an underground detector, because the 
muons they produce will be below the threshold energy. 

There is a slight dip in this distribution as X approaches A'jy^, which results from the approximation 
made to join the high and low energy solutions for the approximation to Eq. \TS\ The low energy 
solution will go to zero when X = and below zero when Eth ^ e^^. The reason for this is 
that these low energy muons have such little energy that they decay in flight, producing a deficit 
in muons anticorrelated to positive temperature changes (the "negative temperature coefficient" 
related in older literature HIIHIISID . This effect is not seen by detectors deeper than 50mwe. The 
fact that there is a dip and subsequent rise in Fig. [T]for X > 480 g/ cm^ does not affect an analysis 
for a detector deeper than 50mwe since the weight is integrated over the entire atmosphere in 
discrete steps of dX and properly normalized, so this atmospheric depth is unimportant for the 
production of relevant muons. 

Remembering that r](X) = (T(X) — Tefr)/Tcfr, the relationship between atmospheric temperature 
fluctuations and intensity variations can be written: 



AI,= A-^dE^= dXaiX)—^, (20) 

J£'th Ct-C/^ Jo Jeff 

where the temperature coefficient a{X) can be written: 

a{X) = {l-X/^^fe-""^- r dE. . . n .2 = W'\X)E;^'^'\ (21) 



with Vr*^(X) given by: 



WiX)^ — ' ^. (22) 

l + {l + l)Bl,K{X){E,^co^eiel,f 

The approximation to the integral follows from arguments made by Barrett |l2l. The derivative of 
this expression agrees with the integrand in the energy region of interest to within 2%. The weights 
as a function of X using the threshold energies of the detectors under consideration can be seen in 
Fig.[il The fact that the lines are nearly on top of each other shows that the weight of the particular 
atmospheric depth does not depend very much on the threshold energy. Recalling that M applies 
equally to K and vr mesons and that the total muon intensity is the sum of the contribution by K 
and TT (Eq.[T2l), the temperature induced change in muon intensity can be written: 



A/, = rdXa-{X)^^ + rdXa^iX)^^. (23) 

Jo leS Jo leS 



Letting T^g be defined such that when T{X) = T^s, AI^ = gives 



/q- dXT{X)a%X) + /q- dXT{X)a^{X) 

J^dXa^{X) + J^dXa^iX) ' ^ ^ 
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Since the temperature is usually measured at discrete levels, the integral is calculated numerically 
over the atmospheric levels AX„: 



Tos ^ . ...rl . ..rr.. ■ (25) 



W^^^ is W^^^ evaluated at X„, 1/A;_^ = l/A^v - 1/K,k, = 120g/cm2, = IGOg/cm^ 
and A/^ = 180g/cm2 [[ni. 



Relative Weight 

0.2 0.4 0.6 0.8 1 




T(K) 

Fig. 1. The average mid-latitude summer temperature at various atmospheric depths (dashed line). The 
vertical range is from 1000 hPa (1 hPa = 1.019 g/cm^), near Earth's surface, to 1 hPa (nearly 50 km), near 
the top of the stratosphere. The solid lines are the weight as a function of atmospheric depth used to find T^q 
(Eq.[22l). The blue lines used the OPERA threshold energy, the red lines used the Super-K threshold energy, 
the black lines used the MINOS FD threshold energy, the green lines used the IceCube threshold energy, the 
violet lines used the Baksan threshold energy and the magenta lines used the MINOS ND threshold energy . 



With this definition of Effective Temperature, an Effective Temperature coefficient, ar, can be 
defined: 

1 r roo poo 

/ rfXa^(X)+ / dXa^{X) , (26) 
Jo Jo J 



"^ = 70 



where is the intensity for a given temperature T. Now that the atmospheric temperature has been 
parametrized and ar defined, the relationship between atmospheric temperature fluctuations and 
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Normalized Intensity 

0.2 0.4 0.6 0.8 1 




T(K) 



Fig. 2. The solid lines show the meson intensity as a function of atmospheric depth (Eq. \^ and the dot- 
dash lines show the muon intensity as a function of atmospheric depth (Eq. IT). The normalization in Eq. |4] 
(NQZ]\fM /^n) was set to 1 to show the dependence on X more clearly. The range of the expressions were 
adjusted so that the maximum value of the OPERA meson intensity was equal to 1, and the other equations 
were scaled appropriately. These figures were produced with particular energy values corresponding to the 
threshold energy of the detectors under consideration. The blue lines used the OPERA threshold energy, the 
red lines used the Super-K threshold energy, the black lines used the MINOS FD threshold energy, the green 
lines used the IceCube threshold energy, the violet lines used the Baksan threshold energy and the magenta 
lines used the MINOS ND threshold energy . 



intensity variations can be written: 



^ . f ,XaiX)^ ^ .A (27, 

J,, -'0 icflf Jeff 



Note that the expression to calculate T^s in the pion scaling limit, ignoring the kaon contribution 
is the same as the MACRO [3] calculation for effective temperature. This distribution reflects the 
dominant atmospheric phenomena that produce muons visible to a detector underground. High 
energy mesons produced at the top of the atmosphere have the greatest influence on the seasonal 
variation because they are created where the density is lower, so they have the highest probability 
to decay into muons. High energy mesons that are produced lower in the atmosphere have a greater 
probability of interacting a second time, and thus greater probability of producing muons that are 
not seen by an underground detector. 
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4 Theoretical Effective Temperature Coefficient 



The theoretical prediction of ax for properly weighted atmospheric temperature distribution can 
be written (Eq.lTTi): 



JO dT ■ 



(28) 



Barrett |l2l shows that for a muon spectrum such as Eq. [8l the theoretical «r can be written: 



aT 



-7. 



(29) 



The prediction for can be calculated using the intensity found in Eq.[T2land a little algebra: 

1 l/e'K + A'ADjDKr/el 



where 
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(31) 



7 + 1 l.l^th cos 61 

Note that this can be reduced to MACRO'S previously published expression (a^)^ which was 
only valid for pion induced muons, by setting A]^ = (no kaon contribution). This approximation 
can be extended to a kaon-only temperature coefficient, {aT)K by setting the pion term (first term 
in parentheses) in Eq. [121 The result is an independent model of the temperature coefficient for 
each of the meson species: 
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7 + 1 l.l^^th cos6' 



+ 1 



(32) 



To compare the experimental ar to the theoretical expectation, a simple numerical integration 
using a Monte Carlo method was performed. A muon energy and cos 9 were chosen out of the 
differential muon intensity (Eq.[T2l). A random azimuthal angle, 0, was chosen and combined with 
cos 6*. Column depth was calculated as d = h/ cos 9, where h is the detector depth in mwe for 
standard rock with flat overburden. The threshold surface energy required for a muon to survive 
this column depth is found from the expression for threshold energy (Eq.fTTI). If the chosen i?^ was 
greater than Eth, it was used in the calculation of the theoretical (ay). This was repeated for 10,000 
successful muons with > Eth, at depths from to 4,000 mwe. The result of this calculation, 
along with the experimental result from the MINOS experiment [HI can be seen in Fig. [3] as the 
solid line. This curve includes the "negative temperature effect" (muon decay correction) term, 
5' = cos ^)(m^c2/7/cr^) (7/7 + 1) ln(1030/A^ cos e) fT\, which goes to zero for E^, > 

50 GeV. 

The kaon component of air showers that can be observed at 1400 mwe is about 10%, but the energy 
is too low for kaon-induced muon production to be affected by changes in temperature. The result 
is the large gap between the pion only curve and the Kvr curve. As the depth increases, the energy 
of sampled muons also increases, which results in a greater contribution to ar by kaon induced 
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% 1000 2000 3000 4000 

depth (mwe) 

Fig. 3. The theoretical {aT{X)) (solid curve), the {ax {X)) ^{dashed curve), the (aT(-'^))x(dotted curve) 
for slant depths up to 4,000 mwe. The MINOS data point is from [5 1. The cross-hatched regions indicate the 
sensitivity (separation between the three models for a particular depth) that current underground detectors 
have to measurements of (axiX)). 

muons. The asymptotic behavior of the theoretical ax approaching one as primary energy increases 
is expected from Eql28l At very high primary energies, the intensity is proportional to the critical 
meson energy, which depends on temperature. Thus, for an isothermal atmosphere, intensity will 
be directly proportional to the temperature (the constant of proportionality, ax, will be one). 



5 Method for Measurement of Atmospheric K/tt Ratio 



The theoretical uncertainty of the atmospheric K/vr ratio is of order 40% [16J. There was not 
a measurement of this ratio with cosmic rays until 2009, when it was made using MINOS FD 
data im. Previous measurements have been made at accelerators for p+p collisions iflTll . Au+Au 
collisions [HI, Pb+Pb collisions [[T9l20ll . 



The expression for the experimental ax is written in Eq. [271 The kaon influence causes an overall 
decrease in total ax, shown in Sec.[4l Because the left hand side of Eq.[27]depends only on counting 
rate, it can be broken into meson components: 

^ ^ - ottt^, (33) 



n 




which can be rewritten: 

Recall that in the pion scaling limit, only pions are assumed to contribute to the seasonal effect. 
From that, a model for pion-only and kaon-only temperature coefficients were developed in Eq.[32l 
Such a seasonal effect can be written: 



The ratio of the muon counting rates / RJ^ is equivalent to the ratio of muons from kaons to 
muons from pions Njf /NJ^, which will be written r^(K/n). Rearranging and inserting Eq.[35]for 
both kaons and pions into Eq. [33] gives: 



r^{K/7r) = — (ar)^ + (aT)i^ " 1 (36) 



_ {aT)w/aT - 1 

1 - {aT)K/0iT 

The value for r^{K/'n) can be predicted by integrating Eq 

J 

where 



(37) 



r^{K/7:) = ^ = C2xA],, (38) 



^ JBthCose 1 + cost^/eii' 

n= / —rdE. cos 9. (40) 

^ cose 1 + I.IE^ cos ^/eO ^ ^ ^ 

The parameter A]^ is defined as 

0.635 x.(A7.)|i^le^, (41) 

where r{K/n) = [fTTIl is the ratio of kaons to pions produced in the primary cosmic ray 
interactions. Inserting Eq.l4T]into Eq.[38]and rearranging gives an expression for r{K/n) in terms 
of r^{K/n): 

r(KM ^ i- X r,(KM x (42) 
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The MINOS-FD (5) measurement is shown in Fig. H, along with STARS [[H, NA49 II19I20II and 
ISR [17J accelerator measurements. The cross-hatched regions show the energy regimes to which 
the underground detectors discussed in this paper are sensitive for K/n ratio measurements using 
the method described. An OPERA measurement would extend to the region an order of magnitude 



^ 0.2- 
0.15 
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▼ STARS (Au+Au, 
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MINOS-ND Baksan IceCube 
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100 



200 
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^(GeV) 



Fig. 4. A compilation of selected measurements of r(K/7r) for various primary particle center of mass ener- 
gies {^/s). The STARS value was from Au-i-Au collisions at RHIC flSl . the NA49 measurement was from 
Pb-i-Pb collisions at SPS II19I20L the ISR measurement was from p-np collisions ITtI, and the MINOS value 
was from cosmic ray primaries + atmospheric nuclei collisions fS^.The thick horizontal bars near the bottom 
of the graph show the typical ranges of cosmic ray primary energies for the collisions that produce muons 
observed by the underground detectors indicated. 



beyond the energy of current fixed target experiments. With a detector area roughly 1,000 times 
the area of the MINOS FD, IceCube should have a cosmic ray muon rate of upon completion of 
construction (Boreal Spring, 2011) of 1,700 Hz [|21||. Assuming a temperature data set of compa- 
rable quality to the BADC ECMWF data (221 used by the MINOS-FD analysis [5], the statistical 
uncertainty in ax could be reduced to ±0.001. Since the absorber material surrounding IceCube 
is ice instead of rock containing iron veins, the column depth should be more well known. This 
could reduce the uncertainty in the depth map, the dominant source of uncertainty in {aT)TT,K, by 
half. These factors taken together could reduce the uncertainty in the measurement of r{K/n) by 
16-30%. 
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6 Summary 



A new method was developed to include the effect of kaons in measurement of seasonal variations 
in underground muon intensity. A temperature coefficient that accounts for the kaon contribution 
was described, and a kaon-inclusive model was defined. These methods were applied to MINOS- 
FD data [5], and the new model fit the data better than the pion only model [3]. A formula was 
described so that other underground experiments, OPERA, Super-K, IceCube, Baksan and the MI- 
NOS ND could quickly apply this method to their data. Pion and kaon decay are affected differently 
by temperature variations, and this difference suggested a method to measure the atmospheric K/vr 
ratio. This method was developed and a formula was offered for other underground experiments to 
follow. 
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